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Diffusion is often accompanied by a reaction or sorption which can induce temperature inhomo- 
geneities. Monte Carlo simulations of Lennard- Jones atoms in zeolite NaCaA are reported with a 
hot zone presumed to be created by a reaction. Our simulations show that localised hot regions can 
alter both the kinetic and transport properties. Further, enhancement of the diffusion constant is 
greater for larger barrier height, a surprising result of considerable significance to many chemical 
and biological processes. We find an unanticipated coupling between reaction and diffusion due to 
the presence of hot zone in addition to that which normally exists via concentration. 

PACS numbers: PACS Numbers: 05.40.Jc,05.60.Cd,82.20.Wt,82.33. Jc,82.75.Mj 



Diffusion within porous materials or confined geometry 
is still poorly understood jj], despite increased attention 
in recent times [FJ, ||. Life sciences has a number of in- 
stances which relate to diffusion within confined regions 
- for example, ion diffusion across membranes and ap- 
proach of a substrate towards an active site of an enzyme 
H . Hydrocarbon separation and catalysis within zeolites 
provide instances of processes in chemistry ||. Prob- 
lems involving fluid flow and excitonic transport through 
porous medium are examples from physical sciences [0. 
The richness of the subject partly arises from the ge- 
ometry of confined systems ( the fractal nature of the 
pores, for instance). Further, while non-uniformity of 
concentration has been dealt with in great detail, that of 
temperature has received little attention. In particular, 
when temperature is inhomogeneous, the very definition 
of diffusion as being an activated process needs a gen- 
eralization. Such non-uniformity in temperature arises 
routinely in biological, chemical and physical systems for 
variety of reasons. Here, we discuss issues relating to the 
possible sources of such hot spots and their influence on 
transport properties in the context of zeolites. 

Zeolites are porous solids with pore sizes comparable 
to molecular dimensions j7j . Due to its rich and diverse 
catalytic as well as molecular sieve properties [|| it has 
attracted much attention. The existence of specific cat- 
alytic and physisorption sites coupled with their poor 
thermal conductivity could lead to local hot regions 
(Typically in lOps, the hot region decays less than a few 
percent). This may affect both kinetic and diffusion prop- 
erties. Such a situation can arise in many biosystems as 
well. For instance, plasma membrane protein-encoding 
m-RNA IST2 is shown to have high asymmetry in con- 
centration between the mother cell and the bud|^]. One 
possible way of maintaining such an asymmetry against 
the concentration gradient is through localized hot or 
cold regions. In spite of the importance of such reaction 
induced hot spot and its influence on the diffusion of the 
species, this problem has not been addressed so far. 

Here, we study the effect of inhomogeneous temper- 
ature presumed to be created by a 'reaction', on the 
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FIG. 1: The effect of a hot zone at BC in the potential ABCD 
is to lower D to D' . 



equilibration rate and self diffusion coefficient of guest 
molecules in zeolites. Monte Carlo simulations on simple 
argon atoms in zeolite A are reported here. Our results 
show that self diffusion coefficient D is increased sub- 
stantially due to the presence of a hot zone. More signif- 
icantly, at a conceptual level our analysis shows that local 
changes in temperature resulting from reactions can in- 
duce additional coupling between reaction and diffusion. 

Landauer||l0||, in a seminal paper, addressed the effect 
of a nonuniform temperature bath on the relative occu- 
pation of competing local energy minima. For the case of 
a bistable potential U(x) (the curve ABCD in Fig. |lj), he 
showed that the presence of a localized heating in a region 
(say BC) lying between the lower energy minimum A and 
the potential barrier maximum C can raise the relative 
population of the higher energy minimum D over that 
given by the Boltzmann factor exp(—AE/kBT). This 

icj. Since 
andauer, 



has come to be known as the 'blowtorch' effect 
this effect is rather counter intuitive, following 
we convey the basic idea. Consider the motion of an 
overdamped particle in this potential (curve ABCD in 
Fig. ^) subject to a uniform temperature To along the 
coordinate. Then, the probability of finding a particle 
at x is P(x) ~ exp(—U{x)/kBT()). Then the probability 
at A is higher than that at D. If the temperature in BC 
is raised to Tb, then P{x) ~ exp{—U{x)/kBTb) in BC, 
which is much smaller than P(x) with T = Tq and hence 
log(P(x)) is flat in BC. Since —logP{x) can be considered 
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FIG. 2: Two cages of zeolite NaCaA. • show reaction sites. 
The saddle point is at the window. A schematic one dimen- 
sional potential along the z-direction is shown below along 
with the induced hot zone. For set A simulations the poten- 
tial at extremes of z goes to infinity (not shown). For set B, 
the potential shown above is periodic. 



as 'potential' U(x)/kBT, raising the temperature to T, 
in BC is equivalent to modifying the 'potential' to flat- 
ter curve BC (Fig. [I]). Since P{x) is unaffected in other 
regions, the curve outside BC will be the same excepting 
that the curve CD would start at C' and end at D' such 
that U(xc) — U(xd) — U(xc) — U(xD')- Thus, the point 
D is brought down relative to A. Consequently, P(xd) 
is higher than at the lower minimum xa- 

A decade later an appropriate diffusion equation for 
a nonuniform temperature profile was derived by van 
Kampen|ll[ , followed by additional work from Landauer 
himself |l2|. Recently, Bekele et aZ.jl3| have shown that 
the escape rate is enhanced due to the presence of a hot 
zone. This implies that the rate of surmounting the bar- 
rier in zeolite will be increased by the presence of a hot 
zone and a consequent increase in diffusion constant. 

The physical system we simulate consists of NaCaA 
zeolite with large (« 11. 5A dia) cages (the supercages) 
interconnected by shared narrow 8-ring windows (« 4.5A 
dia). The potential energy landscape has a maximum 
near the 8-ring window and a minimum located deep 
within the supercage. A species arriving at a hetero- 
geneous reaction site, assumed to be located between the 
window and the center of the cage, (Fig. ||) releases a heat 
q creating a local hot zone. Consequently, the molecule 
surmounts the barrier more easily. Other molecules be- 
hind it also cross the barrier with relative ease due to the 
hot zone. Here, we mimic the reaction by its principal 
effect - the presence of a hot zone - by introducing it in 
between the potential maximum and minimum. In case 
of physisorption we assume that the heat released due to 
sorption is taken up by another diffusing particle. The 
effect arising out of the occupied volume of the adsor- 
bent is ignored as it is negligible (for methane-NaY, this 
is typically 2% of the cage volume) . 

We consider Lennard- Jones particles confined to 2 x 



2x1 unit cells of NaCaA zeolite. There are 2 3 cages 
in each unit cell and 4x4x2 cages in the x, y and z- 
directions. The potential along the z-direction is a sym- 
metric double well(Fig. ^). Since the rate determining 
step is the passage through the 8-ring window, the dis- 
tance of the particle from this plane may be considered as 
the diffusion coordinate. Due to lack of appropriate tech- 
niques for including non-uniformity in temperature (that 
is, in maintaining steady state excess temperature in a 
localised region) within the existing molecular dynamics 
algorithms, we use Metropolis Monte Carlo algorithm in 
the canonical ensemble where the total energy is 

* = ££^££^ a) 

a z gg 

and <t>ab = 4e a6 [(er Q6 /r) 12 - (o- ab /r) 6 ]; a, 6 = g, h is the 
Lennard- Jones potential. Thus, the particle diffuses on 
the potential energy landscape of the zeolite. The poten- 
tial parameters for the guest-guest and guest-host inter- 
actions are a gg — 2.73 A and e gg — 0.9977 kJ/mol. cto-o 
= 2.5447 A, a Na - Na = 3.3694 A, a Ca -Ca = 3.35 A, e -o 
= Z*1.28991 kJ/mol, e Na - Na = ;*0.03924 kJ/mol and 
eca-ca = Z*9.5451 kJ/mol @ where I = 2,4. The cross 
terms are obtained from Lorentz-Berthelot combination 
rule. The cut-off distance was 12 A. 

Two different sets of simulations A and B are carried 
out. In set A, for the calculation of the escape rate, we 
impose periodic boundary conditions(PBC) along the x- 
and y-directions and a repulsive potential (1/r 12 ) at both 
ends along z-direction which enables a comparison with 
earlier workjl3|]. In set B, for the calculation of D in 
three dimensions, PBC along the all the three directions 
is essential. We investigate the influence of degree of 
hotness defined by s = (TJ,— Tq)/Tq and the barrier height 
U a (cx le) on the equilibration rate and D. Here, To 
and Tb are the background and blow torch temperatures 
respectively. To is kept at 300 K, and T, is varied. 

Initially, all the 64 guest particles, at a concentration of 
2 atoms/supercage are uniformly distributed in the four 
left cages located along the z direction (see Fig. |l|). It 
is well known that transport properties are determined 
by the time scales associated with the approach to the 
steady state. These are obtained by allowing the sys- 
tem to evolve towards the steady state in the presence 
of the hot zone. Let ni and n r be respectively the total 
number of particles in the cages to the left and the right 
of saddle point along the z-direction. Then, the decay 
rate is obtained by plotting the fraction of particles in 
the left cages as a function of time. A typical plot of 
In (ni/(ni + n r )) vs. t is shown in Fig. |[ The slope then 
gives the equilibration rate. The curve shown in Fig. |3| 
will reach a plateau at much longer time scale. 

We consider the influence of the hot zone on the rate 
of approach to the steady state as a function of s when 
the hot zone of width w = 2 A is placed at a distance 
d = 1.2775 A from the window. Instead of the escape 
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FIG. 3: A typical plot of ln(ni/(ni +n r )) vs. MC steps 
obtained from MC runs(set A). 
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FIG. 4: A plot of the enhancement factor, /i, as a function 
of the degree of hotness, s = (Ti, — To) /To, for two different 
values of well depth, le (I = 2(D), 4(o)) from MC runs(set A). 



TABLE I: Dh/Do for different sets of parameters. Here d = 
1.2775 A and w = 2 A obtained from MC runs(set B). 
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1.5088 


1.4283 


3 
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3.5416 


2.2703 


2.1971 
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2.6132 


1.8066 


1.6747 



rate, we use the enhancement factor fb, which is the ra- 
tio of the escape rates from the left cage with and with- 
out the hot zone. A plot of as a function of s, for 
two values of I is shown in Fig. || The tendency to ap- 
proach saturation is seen for both values of U a (oc le), 
even though, it is less pronounced for the higher value 
(/ = 4). More importantly, f b is a sensitive function of 
the barrier height 1 13]] (Fig. Q). Thus, the enhancement in 
fb is greater when U a is larger which implies larger en- 
hancement in D for systems with higher energy barriers. 

The significance of these results becomes apparent on 
examining a real system such as methane in faujasite. 
The energy at the physisorption site for methane in NaY 
zeolite (with Si/Al = 3.0) is -18 kJ/mol, the energy differ- 
ence between a free methane and a physisorbed methane. 
However, the energy released is significantly lower( ~ -6 
kJ/mol) when it is already within the zeolite [fPo]. This 
can raise the temperature in the vicinity of the site. 

We now consider the influence of a hot zone on dif- 
fusion through set B simulations for which the starting 
configuration is the final configuration of set A. We have 
studied the influence of s, and the barrier height. The 
ratio of the diffusion constant, Dh with the hot spot to 
that without, D is enhanced in each case(Table ||). We 
note that since all other conditions of the simulation are 
identical in the two situations, D^/Dq is independent of 
the basic length and time scales and also the details of 
the simulation such as the particle displacement (similarly 
for set A.) Note that the larger the barrier height greater 
is the enhancement in D, a result that has considerable 
implication. These results are better understood on the 
basis of simple arguments to estimate Dh/D . 

In the Kramers picture, the slowest time scale deter- 
mining the approach to the steady state is identified with 
the escape rate. To estimate this, consider a one di- 



mensional symmetric double well potential with PBC. 
The rate equation for the number of particles to the 
left (ni) and right (n r ) of the potential maximum, is 
hi = —wi r ni + w r in r , where wi r and w r i are the escape 
rates from the left to the right well and vice versa re- 
spectively. These are also the very time scales that de- 
termine the steady state through detailed balance condi- 
tion: hi = or wi r ni = w r in r . In the absence of hot spot, 
we have wi r = w r i. Then the mean escape rate from a 
well is wo = {wir + w r i)/2. Then the diffusion constant 
in one dimension, Dq = a 2 u>o/2, where a is the distance 
between the two minima. On introduction of a hot zone 
in the left well, wi r has been shown to increase consid- 
erably, while w r i increases only marginally |T^| , and to a 
good approximation w r i ~ wq- The mean escape rate is 
Wh = (wir + w r i)/2. Using this with u>i r — w r in r /ni, we 
get Dh = a 2 Wh/2 ~ D (l + n r /n{)/2. Since, n r /ni > 1 
in the presence of a hot zone, it is evident that D is 
enhanced. Note that this refers to a nonequilibrium in- 
homogeneous situation as n r /ni refers to the steady state 
which can only be obtained numerically. Table |, shows 
the values of Dh/Do obtained from the mean square 
displacement through MC simulations for three sets of 
parameter values along with Dh/Do estimated from the 
above expression using the steady state values of n r /ni 
from MC simulations. (Note n r and n; refers to the total 
number of particles to the left and right of the potential 
maxima between the two cages in the z-direction.) Note 
that the estimated values are close to MC values. 

To facilitate comparison of these results with real sys- 
tems, we have estimated the likely increase in temper- 
ature when hydrocarbons and other guest species are 
sorbed within zeolites such as NaX. We have listed in 
Table || isosteric heat of sorption (AH a d s ) of some linear 
alkanes, Xe and water within faujasites. We have also 
listed the mean heat capacities (C m ) of the guest-zeolite 
systems 0. From these data, the maximum increase in 
AT can be estimated from AT = T b - T = AH ads /C m 
which is in the range 820 K to 2028 K (Table ||). Thus, 
the parameter s varies from 1.7 to 6.7 for which fb can be 
as large as 3 for s ~ 6 and 1 — 2. Since fb is determined 
by the very rate constants that lead to steady state, it 
also implies that the diffusion constant in the inhomoge- 
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TABLE II: Expected rise in temperature and s for typical 
guests when adsorbed in common zeolites estimated from heat 
of adsorption, AH a ds and the mean heat capacity, C m data. 



System 


A U 

L -^ J - J -vap 

kj/mol 


kj/mol 


J/mol.K 


J-b J-0 

K 


s 


guest 


zeolite 


n-C4Hi 


Na-X 


66 


174 


105 


1689" 


5.6 


n-C 7 Hi 6 


Na-X 


87 


228 


176 


1809 6 


6.0 


n-C 7 Hi 6 


Na-X 


87 


228 


209 


1090 c 


3.3 


neo-C 5 Hi2 


Na-X 


54 


130 


129 


1011 6 


3.3 


iso-CgHig 


Na-X 


88 


246 


185 


1329 d 


4.0 


Xe 


Na-Y 




18 


22 


820 e 


1.7 


H 2 


Na-X 




142 


70 


2028 6 


6.7 



"Calculated from AH vap and the ratio of AH a ds to 
AH vap $. b T =300K; c T =333K; d T =325K; e T =473K. 



neous medium can increase by a factor of two, even for 
moderate values of s as can be seen from Tabic |. 

An interplay of reaction and diffusion is known to give 
rise to complex dynamics which can manifest in different 
ways [fl6| . The product profile in a reaction is controlled 
by the diffusion rate of the product species formed rather 
than the reaction rate as is the case with the formation 
of p-xylene in ZSM-5Q. Due to low diffusivity of ortho 
and meta isomers, they are not observed as products even 
though they are formed. The coupling between reaction 
and diffusion in such systems is via the concentration of 
the reactant and product species jL7|. In contrast, in 
the present situation, the enhancement in D is a direct 
consequence of inhomogeneous temperature. This study 
demonstrates that such a coupling between reaction and 
diffusion can arise not just due to concentration, but also 
due to the increase in local temperature, a fact that could 



not have been anticipated. These results also show how 
diffusion is increased in the presence of physisorption or 
chemisorption which are usually exothermic. D may de- 
crease if the reaction is endothermic. 

The present analysis can provide an insight into a well 
known experimental observation where a warm adsorp- 
tion front is seen to move rapidly during the adsorption 
of a gas into an evacuated zeolite Jl8[. As the initial 
molecules arrive at a physisorption site, heat is released 
which aids the molecules at the front to cross over the 
energy barrier and propel the gas forward. Zeolites are 
crystalline solids, and hence the active sites are located 
in a periodic manner. Thus, as the front moves further 
into the zeolite, hot zones are created successively pro- 
viding a periodic driving force for diffusion which is in 
addition to that arising from the concentration gradient. 
This explains the rapid movement of the front. 

Thus, influence of a reaction induced hot zone can af- 
fect diffusion in ways which will be important both from 
fundamental as well as industrial perspective. Though 
our results arc obtained in the context of zeolites, it is 
evident that they are of significance to many biological 
processes where concentration gradients are frequently 
accompanied by difference in temperature. These results 
also have implication to the separation of isotopes fll9|, 
the petrochemical industries, fast ion conducting battery 
materials, etc. We believe that in the foreseeable future, 
local inhomogeneities in temperature will be exploited in 
a number of ways to bring forth novel processes. 
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